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Key Points: 


¢ Joint inversion of Swarm and CHAMP data from magnetospheric and ocean tidal cur- 
rents reveals global features of the upper mantle and MTZ. 

- Integrated imaging of the mantle electrical conductivity from multiple sources signifi- 
cantly improves resolution. 

¢ Retrieved conductivity profiles are compatible with a pyrolytic composition and a 
moderate mantle water content. 


Corresponding author: A. V. Grayver, agrayver@erdw.ethz.ch 


27 


28 


29 


30 


31 


32 


33 


34 


35 


36 


37 


38 


39 


40 


41 


42 


43 


44 


45 


46 


47 


48 


49 


50 


51 


52 


53 


54 


55 


56 


57 


58 


59 


60 


61 


Abstract 

We present a new global electrical conductivity model of Earth’s mantle. The model was 
derived by using a novel methodology, which is based on inverting satellite magnetic field 
measurements from different sources simultaneously. Specifically, we estimated responses 

of magnetospheric origin and ocean tidal magnetic signals from the most recent Swarm and 
CHAMP data. The challenging task of properly accounting for the ocean effect in the data 
was addressed through full three-dimensional solution of Maxwell’s equations. We show that 
simultaneous inversion of magnetospheric and tidal magnetic signals results in a model with 
much improved resolution. Comparison with laboratory-based conductivity profiles shows 
that obtained models are compatible with a pyrolytic composition and a water content of 0.01 
wt% and 0.1 wt% in the upper mantle and transition zone, respectively. 


1 Introduction 


Electromagnetic (EM) sounding is an important technique for studying Earth’s interior 
and its material properties. It can be used to infer electrical conductivity in depth and map 
its lateral variations within the Earth, thereby carrying information about composition, tem- 
perature and presence of water or melt in the mantle [Karato, 2011; Katsura and Yoshino, 
2015; Khan, 2016]. The unique characteristic of EM methods is the wide frequency range 
, corresponding to sounding depths from crust to lower mantle. However, across this wide 
frequency range several excitation mechanisms co-exist [Kuvshinov, 2008]. This requires ad- 
justments in the source parametrization during data processing and modeling stages. There- 
fore, interpretation of these data is usually done separately, resulting in reduced resolution of 
individual models and potentially causing inconsistencies between them. 


Simultaneous inversion of multiple data sets from different sources should lead to im- 
proved resolution and smaller uncertainties, and, as a consequence, additional constraints on 
the fundamental aspects of the composition, structure, and dynamics of the Earth. For in- 
stance, Egbert et al. [1992] and Bahr et al. [1993] estimated responses using ionospheric and 
magnetospheric signals, and obtained regional conductivity models of the upper and lower 
mantle with better resolutions due to wider frequency range of the combined responses. 
However, these studies used only land observatory data. Since then, operation of low-orbit 
satellites (Oersted, CHAMP, SAC-C, Swarm) [Olsen et al., 2013] has provided a wealth 
of data opening new opportunities for mantle conductivity studies. In contrast to land ob- 
servatories, processing of data coming from constantly moving satellites is more challeng- 
ing, since it requires a careful separation of the magnetic fields from different sources. Re- 
cently, new data processing approaches have experienced a significant progress [Sabaka 
et al., 2015], enabling the usage of both time-varying magnetospheric and tidal magnetic 
fields for global EM sounding. However, sensitivity of the methods depends on frequency 
content and mechanism of excitation. For instance, long period (periods > 1.5 days) mag- 
netospheric responses are more sensitive to the conductivity in the mantle transition zone 
(MTZ) and below [Kuvshinov and Olsen, 2006; Velimsky et al., 2006; Civet et al., 2015], 
whereas tidal magnetic signals , specifically signals due to lunar principal semi-diurnal M> 
tide, are more sensitive to upper mantle conductivity [Grayver et al., 2016]. As a conse- 
quence, simultaneous inversion is expected to provide improved depth resolution. To cor- 
roborate this here, we perform inversions of magnetic data derived from the most recent 
satellite data and compare individual and joint inversion results, in addition to comparison 
with laboratory-based conductivity profiles for the purpose of making thermo-chemical infer- 
ences. 
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2 Methods 
2.1 Satellite data 
2.1.1 Magnetospheric responses 


For periods longer than one day, signals due to magnetospheric ring current dominates 
the measured time-varying magnetic fields [cf. Piithe et al., 2015a]. These signals are con- 
ventionally described by the first zonal spherical harmonic. In this work, we derived magne- 
tospheric responses through the so called QO’ (w)-response [e.g. Piithe and Kuvshinov, 2013), 
which relates frequency-dependent inducing, Ee) (w), and induced, i’(w), coefficients as 


iw) = O°(w)e(w). 


From this, the global C,-response [e.g. Olsen, 1999] on the surface of the Earth can be 
calculated as 


gle 208 (w) 


C = -——_——_: 
(wv) = 5 1+ 0%) 


Note that for a radially homogeneous Earth, C)-responses exhibit monotonic growth 
with respect to period (T = an). 


To quantify the degree of correlation between the inducing and induced coefficients, 
we used squared coherence given by 


(77), fF (@)) I? 


2 
coh*(w) = ————_____—_—_—__., 
(i)(), 77 @)){ €?(@), 6? @)) 
where (, ) stands for inner product between two vectors. In this context, vectors are given by 
a set of the Fourier-transformed windows of V(t) and Ee (t) time-series. The closer this value 
to its upper bound of one, the more variability in i?(w) can be explained by the variability in 


Ej (w) 


2.1.2 Tidal magnetic signals 


The tidally-induced flow of the electrically conductive ocean water in Earth’s main 
magnetic field generates electric currents, which in turn induce secondary EM field in the 
subsurface as a result of which the total magnetic field measured on land or at a satellite car- 
ries information about subsurface electrical structure. In contrast to other conventional EM 
sources of ionospheric and magnetospheric origin, which are inductively coupled with the 
Earth, the unique characteristic of the motionally-induced ocean currents is its galvanic cou- 
pling with the Earth. This enhances sensitivity to the resistive subsurface structures since the 
induced fields are influenced by the toroidal (galvanic) part of the tidal primary EM field. 


Despite small amplitude, tidal magnetic signals due to the semi-diurnal lunar M) tide 
(period of 12 hours and 25 minutes) have been reliably extracted from satellite measure- 
ments using the Comprehensive Inversion approach based on the simultaneous robust least- 
squares estimators of different contributions (core, crust, etc.) and careful pre-selection of 
data [Sabaka et al., 2015, 2016]. They were used to retrieve upper mantle conductivity under 
the oceanic crust [Grayver et al., 2016]. The inverted signals are represented by the radial 
magnetic field component By? at the satellite altitude. 
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Figure 1. Model parametrization adopted in this study. The model consists of a laterally-varying top-most 


conductivity layer and a number of laterally-homogeneous conductivity layers underneath. 


2.2 Forward modeling 


In this work, we focus on determining the radial conductivity structure under the oceans 
and continents. However, to accurately calculate electromagnetic responses due to magneto- 
spheric or tidally-induced oceanic currents, it is essential to account for non-uniform oceans 
[Everett et al., 2003; Kuvshinov, 2008]. To this end, we added a heterogeneous conductiv- 
ity layer corresponding to oceans and continents on top of the laterally homogeneous model 
(Figure 1). Calculating EM field for such a 3D model requires solution of Maxwell’s equa- 
tions 


Uy V x B = cE + 7 


: 4 
VXE = iwB 4) 


where E and B are electric and magnetic fields, respectively; jo is magnetic permeability of 
vacuum; o electrical conductivity; w the angular frequency and j**' the extraneous current. 
We assume e~'”! sign convention. 


To solve system (4) numerically, we used global solver [Kuvshinov, 2008] based on the 
integral equation approach. 


For tidal flow, the extraneous current is confined to the oceans and is given by 
F*(, 8) = o5(4, 9) (8(4, 8) x B™(4, 6)), (5) 


where o-, is the conductivity of seawater, Bin is Earth’s main (core) magnetic field, v = 

u/h, h is the height of the water column and u is the depth-integrated seawater velocity due 

to tidal forces. Symbols ¢ and 6 denote , respectively longitude and co-latitude. See Grayver 
et al. [2016] for more details about eq. 5 individual terms. 


For the global QV (w) response, which we need to derive the global C\-response, the 
extraneous source current is parameterized using a single S°(@) = cos @ spherical harmonic. 
The source is then represented as a current sheet located above the Earth’s surface. Once sys- 
tem (4) is solved for the given current distribution, and the radial component of the magnetic 
field, B,, at the Earth’s surface is obtained, the O°(w) is expressed via surface integral in 
geomagnetic coordinates as 


3 x. = 4 
OH) =F ff (Bw. — Be" w.F)) SP Gores, 6) 
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where B&*' is the external magnetic field, 7 = (r = a,¢Gm,.9Gm) is the position vector in 
geomagnetic coordinates on the surface of the Earth, respectively, and a = 6371.2 km is the 
mean radius of the Earth. 


2.3 Stochastic inversion of multi-source data 


The unknown conductivity values a --- oy (Figure 1) can be estimated from satellite 
responses by solving a non-linear inverse problem, which we formulate as a minimization 
task 


M 
_ {[¢atm) £6 
aremin | ——— + — ];m|?” }, 
gn ; bm a | 


where m = [A(o1)---A(oy)] € R” is the vector of unknown model parameters and A(-) 
represents a log-based transformation ensuring positivity of the argument [e.g. Key, 2016]; 

f is a regularization parameter; I; is a regularization operator for the i-th model parameter; 
and scalar p,, controls the norm of the regularization term. By varying py», one retrieves 
different regularization norms, ranging from smooth Lyz-norm (p,, = 2) to structurally sparse 
L\-norm (pm = 1) solutions. Special attention is paid to the data misfit term given by 


Nx 


baim) = Y7 | — wk (fam) - AP), 
1 


kEM k i= 


where M is a set of methods and w*, f*(m), and d* are corresponding data weights (recipro- 
cal of uncertainties), forward operator, and observed data, respectively. Note that normaliz- 
ing with the number of actual measurements (;,) is an important aspect that helps balance 
contributions of different methods in the total misfit term of the minimized functional. In 
general, the approach can be extended to any number of methods, but here is limited to meth- 
ods discussed in Section 2.1. 


Finally, the minimization problem (7) is solved by using a stochastic optimization algo- 
rithm as described in Grayver and Kuvshinov [2016]. 


3 Results 
3.1 Satellite data 


To estimate global C}-responses, we used satellite magnetic measurements. The re- 
sponses were derived from 37 months (from Dec 2013 to Jan 2017) of Swarm data for peri- 
ods of 1.5 - 87 days (Figure 2). For periods > 90 days, we took responses derived from the 
much longer CM5 (combined CHAMP, Oersted and SAC-C data) time-series [Sabaka et al., 
2015]. In order to better account for the complexity of the source, the magnetospheric time 
series were parametrized using spherical harmonics up to degree n = 2 and order m = 1, al- 
though only the term corresponding to the n = 1, m = O was used to estimate C}-responses in 
the frequency domain. This choice is justified since this term is dominant [e.g. Shore et al., 
2016] and most sensitive to the radial structure of the Earth [Kuvshinov, 2008], which we 
aim to recover in this study. Figure 2 shows statistically estimated responses, their uncertain- 
ties and squared coherencies. Clearly, using Swarm data results in higher coherency for peri- 
ods up to = 90 days. For longer periods, coherency drops because of still insufficient length 
of the Swarm time series. In contrast, responses estimated from the CM5 data exhibit lower 
coherencies for periods < 90 days, but due to longer time series (~ 12 years), longer periods 
up to 177 days are better resolved . This motivated our decision to combine responses from 
different missions. Additionally, we used magnetic signals due to the semi-diurnal M> lunar 
tide extracted from 12 years of satellite data [Sabaka et al., 2015]. The radial magnetic field 
component (Figure 3) of this signal was used in the inversion. 
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Figure 2. Cy, responses and their squared coherences estimated from Swarm (Nov 2013 - Dec 2016) and 
CM5 data (CHAMP, Oersted, SAC-C). Positive and negative values represent real and imaginary parts of the 


response, respectively. 
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Figure 3. Real (a) and imaginary (b) parts of the radial magnetic field component due to semi-diurnal M> 


tide at 430 km altitude. Standard deviation of the signals is shown in (c), note different scale. 


3.2 Inversion 


In this study, the subsurface was parametrized using 45 layers ranging in thickness 
from 9 km right under the oceans and continents to 120 km at the core-mantle boundary 
where a metal conductor (a = 10° S/m) is assumed. The starting model was a homogeneous 
spherical shell of 0.2 S/m. 


Figure 4 shows models obtained by inverting satellite magnetospheric and ocean tidal 
signals separately and jointly. Notably, inversion of C)-responses fails to recover a prominent 
boundary between the lithosphere and astenosphere, which results from the lack of resolution 
in the upper mantle [Puithe et al., 2015b]. This is not surprising given that the shortest period 
for C;-responses is 1.5 days (Figure 2). In contrast, the conductivity model obtained by in- 
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Figure 4. Global conductivity models derived from separate and joint inversions of satellite data. The C- 
response profile denotes the model obtained by inverting magnetospheric ring current responses and the M2 
model denotes the global profile derived from the magnetic tidal signals due to semi-diurnal M> tide. Joint in- 
versions were performed using smoothing and structurally sparse (Lj -norm) regularization. Individual models 
were calculated with smoothing regularization. For reference, values for the average lithosphere-astenosphere 


boundary under the oceans and mantle transition zone are plotted as dashed horizontal lines. 


verting tidal magnetic signals displays a sharp conductivity increase around the lithosphere- 
astenosphere boundary (LAB) at the depth of 70-80 km, but does not show any large varia- 
tions below ~ 300 km, where it attains a value close to the initial conductivity model. The 
models obtained from the joint inversion of magnetospheric C)-responses and tidal mag- 
netic signals managed to resolve the LAB and at the same time constrain conductivity of the 
mantle transition zone (MTZ) and below. We used different types of regularization norms 
to produce smooth and structurally sparse models . Both models fit data virtually equally 
well, attesting to the non-uniqueness of the inverse problem and data uncertainties. 


Let us now examine the data responses _ these models produce. Figure 5(a) shows ob- 
served C,-responses as well as responses calculated using the models from Figure 4. One 
sees that the responses calculated for the models derived from the inversion of Cj -responses 
alone and the joint inversion model fit data within uncertainties, whereas the M2 model pro- 
duces substantially different responses. . While the real part of C,-responses for the M> 
model is close to the observed data for periods < 10 days, the imaginary part differs for all 
periods. This behaviour is confirmed through synthetic tests (see supplementary material) 
and is to be expected since the Mz model is not forced to fit C}-responses. Further, Figure 
5(b-d) shows absolute residuals between observed and predicted tidal magnetic signals. Here 
we see that the residuals are systematically larger for the C,-response model (Figure 5c), with 
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Figure 5. (a) Observed and calculated global C)-responses for models shown in Figure 4. Positive and 
negative values represent real and imaginary parts of the response, respectively. (b-d) Magnitude of the radial 
magnetic field component residuals between observed tidal signals and their predicted counterparts for the 
models shown in Figure 4: models obtained by inverting tidal magnetic signals only (b), magnetospheric 


C,-responses only (c) and both simultaneously (d). 


differences reaching up 40% of the original signal amplitude. For instance, the residuals are 
large in regions around South Africa, west of Australia, around New Zealand, west of Cal- 
ifornia, south of Alaska. This suggests that the increase in conductivity at the LAB that is 
missing in this model is required to explain the data. Indeed, and as expected, both the M 
and joint inversion models explain tidal magnetic signals equally well (cf. Figure 5b and 5d). 
Note that since joint smooth and sparse models produce virtually identical responses, only 
smooth model responses are shown in Figure 5. 


3.3, Comparison with laboratory-based conductivity profiles 


Joint inversion models seem to constrain upper and mid-mantle conductivities bet- 
ter than individual inversions. Therefore, it is instructive to interpret these models. To this 
end, we compute laboratory-based bulk electrical conductivity profiles using the approach of 
Khan [2016]. Bulk electrical conductivity is estimated from the mineralogy and databases 
of laboratory mineral conductivity measurements. Whereas equilibrium rock mineralogy, 
including elastic moduli and density, is computed by free-energy minimization [Connolly, 
2009] as a function of pressure, temperature, and bulk composition using the thermody- 
namic formulation and data compiled by Stixrude and Lithgow-Bertelloni [2011]. We model 
mantle composition using the NazO-CaO-FeO-MgO-Al1,03-SiO. chemical system; bulk 
rock conductivity and elastic properties are estimated by employing appropriate averag- 
ing techniques. The pressure profile is obtained by integrating the load from the surface. 
We compute bulk electrical conductivity profiles fora pyrolytic mantle and a standard 
temperature of 1390°C at the base of a 80 km thick lithosphere [Katsura et al., 2010]. The 
sublithospheric mantle adiabat is defined by the entropy of the lithology at the base of the 
lithosphere, whereas in the lithosphere, temperature is computed by a linear geothermal gra- 
dient (see supplementary material). Elastic properties and density produced by this thermo- 


231 


232 


233 


234 


235 


236 


237 


238 


239 


240 


241 


242 


243 


244 


245 


246 


247 


248 


249 


250 


251 


252 


253 


254 


255 


256 


257 


258 


259 


260 


chemical model agree remarkably well with PREM (see supplementary material) of Dziewon- 
ski and Anderson [1981]. 
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Figure 6. (a) Comparison of conductivity models obtained from inversion of satellite data and several 
laboratory-based conductivity profiles calculated following the approach of Khan [2016]. The laboratory- 
based profiles are based on the assumption of a pyrolytic mantle and different water contents in olivine (Ol), 
wadsleyite (Wad) and garnet (Gr). (b) Global Cj -responses calculated for the laboratory-based conductivity 
profiles shown in (a). Observed responses are shown with circles. Positive and negative values represent real 


and imaginary parts of the response, respectively. 


Figure 6 shows a number of laboratory-based conductivity profiles calculated for dif- 
ferent mantle mineral water contents and plotted together with the joint inversion results. For 
present purposes, we varied the water content of olivine, garnet, and wadsleyite. The wa- 
ter contents of clinopyroxene, orthopyroxene, and ringwoodite are estimated using the water 
partition coefficients described in Khan [2016], which are based on the measurements of Jn- 
oue et al. [2010] and Férot and Bolfan-Casanova [2012]. As is evident from the figure, a 
dry mantle produces conductivities which are much lower than the conductivity of the mod- 
els obtained from the joint inversion. Moderate amounts of water [Karato, 2011; Khan and 
Shankland, 2012], 0.01 wt% in olivine and 0.1 wt% in wadsleyite, in the upper mantle and 
transition zone results in conductivities which are much closer to the inverted models. An 
increase of 0.01 wt% in the water content of garnet results in higher conductivities through- 
out the upper mantle and MTZ improving the match to the smooth model and observations 
(Figure 6b). However, these differences are likely within the uncertainty of our models and 
should be considered with caution. The conductivity of the lower mantle in the inverted 
models is close to the laboratory predictions. 


While this interpretation is qualitative and a direct inversion in terms of thermo-chemical 
parameters is more appropriate [Khan, 2016], these results stress that conductivity models 
obtained from joint inversion of data from very different sources produce self-consistent 
models. The thermo-chemical modeling combined with laboratory measurements of the 
electrical conductivity further confirms that these models are consistent with plausible man- 
tle properties and moderate water contents, in addition to radial seismic reference models 
(see supplementary material). 
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4 Conclusions 


The inversion of natural source EM data for Earth’s mantle electrical conductivity 
usually relies on single-source data and therefore faces the problem of limited resolution at 
different depths due to limits in frequency range imposed by varying source morphology. 

We showed that inverting data from magnetospheric and ocean tidal sources simultaneously 
yields a consistent conductivity profile of the upper mantle and transition zone. The obtained 
global profile is capable of fitting individual data types as well as separate inversions and 
efficiently exploits sensitivity overlap between different sources. 


The new conductivity profile provides additional constraints on estimations of geo- 
physically relevant mantle properties through comparisons with laboratory-based conduc- 
tivity profiles. Specifically, assuming a pyrolytic mantle composition and the temperature of 
T = 1380°C at LAB we found that a moderate amount of water is necessary to explain the 
observed conductivity values in the astenosphere and MTZ. However, for the upper mantle, 
this profile is more representative of the mantle under the oceans since tidal signals are neg- 
ligible above continents. Taking these points into account, the new model can serve as a new 
reference for studies, which need to account for mantle conductivity such as in space weather 
or oceanography. Finally, the approach of jointly inverting multi-source data can significantly 
help studies that aim at mapping lateral variations in mantle conductivity. 
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